
# 
# This script reproduces and applies code from the replication files of Yimeng Li (2019)
# "Relaxing the No Liars Assumption in List Experiment Analyses"
# DOI: https://doi.org/10.1017/pan.2019.7

remove(list=ls())

# Load libraries:
library(bindata) # required by sim_bounds.R
library(dplyr)
library(forcats)
library(haven)
library(ggplot2)
library(list)
library(here)
library(nleqslv) # required by list_relaxed_liars.R
library(purrr)
library(tidyr)

source("Li_2019_NoLiars/list_relaxed_liars.R")
source("Li_2019_NoLiars/list_types_plot_NL.R")
source("Li_2019_NoLiars/list_create_df.R")


list_dataNL <- read_dta("data/listNL.dta")
list_dataNL <-as.data.frame(list_dataNL)
list_dataNL<- list_dataNL  %>% 
  mutate(treatfac=as.factor(treated),
         IDvar=as.factor(IDvar),
         round.fac=as.factor(round)) |> 
  mutate(value= value-1)

options(digits = 3)

dist_control_Example <- c(303, 679, 120, 40)
dist_treated_X_Example <- c(157, 342, 549, 77, 39)

DiM_X_Example <- ictreg(y ~ 1, data = create_list_df(dist_control_Example, dist_treated_X_Example, J = 3), treat = "treat", J = 3, method = "lm")

unname(DiM_X_Example$par.treat)
unname(c(DiM_X_Example$par.treat+qnorm(0.025)*DiM_X_Example$se.treat, DiM_X_Example$par.treat+qnorm(0.975)*DiM_X_Example$se.treat))

bounds_X_Example <- list_relaxed_liars(dist_control_Example, dist_treated_X_Example, J = 3, affirmative = TRUE)

bounds_X_Example$bounds
bounds_X_Example$CI

FigureNL<- list_types_NL_plot(bounds_X_Example) 

F1<- FigureNL$plot_two_types
F2<- FigureNL$plot_three_types

F1+F2
ggsave("FigureA9.png",
       path = here("figures"),
       dpi=600)